//////////////////////////////////////////////////////////////////
// (c) Copyright 2010-  by Jeongnim Kim
//////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////
//   National Center for Supercomputing Applications &
//   Materials Computation Center
//   University of Illinois, Urbana-Champaign
//   Urbana, IL 61801
//   e-mail: jnkim@ncsa.uiuc.edu
//
// Supported by 
//   National Center for Supercomputing Applications, UIUC
//   Materials Computation Center, UIUC
//////////////////////////////////////////////////////////////////
// -*- C++ -*-
#ifndef QMCPLUSPLUS_TRANSPOSE2D_OPS
#define QMCPLUSPLUS_TRANSPOSE2D_OPS
#include <benchmark/transpose.h>
namespace qmcplusplus
{
  /** transpose functor using a simple loop
   */
  template<typename MT>
    struct TransposeOMP
    {
      inline void operator()(const vector<MT*>& in, MT& out, int ip)
      {
        int npy=out.rows();
        int nb=in.size();
        int i_offset=ip*npy;
        for(int jp=0; jp<nb;++jp)
        {
          int npx=in[jp]->rows();
          int j_offset=jp*npx;
          trans_b(npx, in[jp]->cols(), in[jp]->data() , npy, out.cols(), out.data() , &i_offset, &j_offset);
        }
      }

    };

  /** transpose functor using [ip,ip+n) with a circular indexing
   */
  template<typename MT>
    struct TransposeOMP2
    {
      inline void operator()(const vector<MT*>& in, MT& out, int ip)
      {
        int npy=out.rows();
        int i_offset=ip*npy;
        int nb=in.size();
        for(int j=ip; j<ip+nb; ++j)
        {
          int jp=(j<nb)?j:j-nb;
          int npx=in[jp]->rows();
          int j_offset=jp*npx;
          trans_b(npx, in[jp]->cols(), in[jp]->data() , npy, out.cols(), out.data() , &i_offset, &j_offset);
        }
      }
    };

  template<typename MT>
    struct TransposeOMP_M1
    {

      inline void operator()(const vector<MT*>& in, MT& out, int m, int ip)
      {
        int nb=in.size();
        int npy_b=out.size(1);
        int npx=out.size(2);
        int nxy_in=npx*npy_b;
        for(int k=0; k<m; ++k)
        {
          int i_offset=k*nxy_in+ip*npy_b;
          for(int jp=0; jp<nb;++jp)
          {
            int npx_b=in[jp]->size(1);
            int npy=in[jp]->size(2);
            int nxy_out=npx_b*npy;
            int j_offset=k*nxy_out+jp*npx_b;
            trans_b(npx_b, npy, in[jp]->data(), npy_b, npx, out.data() , &i_offset, &j_offset);
          }
        }
      }
    };

  template<typename MT>
    struct TransposeOMP_M2
    {

      inline void operator()(const vector<MT*>& in, MT& out, int m, int ip)
      {
        int nb=in.size();
        int npy_b=out.size(1);
        int npx=out.size(2);
        int nxy_in=npx*npy_b;
        for(int k=0; k<m; ++k)
        {
          int i_offset=k*nxy_in+ip*npy_b;
          for(int j=0; j<nb;++j)
          {
            int jp=(j+nb)%nb;
            int npx_b=in[jp]->size(1);
            int npy=in[jp]->size(2);
            int nxy_out=npx_b*npy;
            int j_offset=k*nxy_out+jp*npx_b;
            trans_b(npx_b, npy, in[jp]->data(), npy_b, npx, out.data() , &i_offset, &j_offset);
          }
        }
      }
    };
}
#endif
/***************************************************************************
 * $RCSfile$   $Author: jmcminis $
 * $Revision: 4631 $   $Date: 2010-02-18 14:08:30 -0600 (Thu, 18 Feb 2010) $
 * $Id: transpose.h 4631 2010-02-18 20:08:30Z jmcminis $ 
 ***************************************************************************/
